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' High temperature approximation (HTA) is used to describe the phase behavior of polydisperse 

multi- Yukawa hard-sphere chain fluid mixtures with chain length polydispersity. It is demonstrated 
that in the frames of the HTA the model belongs to the class of "truncatable free energy models" , 
i.e. the models with thermodynamical properties (Helmholtz free energy, chemical potential and 
pressure) defined by the finite number of generalized moments. Using this property we were able to 
calculate the complete phase diagram (i.e., cloud and shadow curves as well as binodals) and chain 
length distribution functions of the coexisting phases. 
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I. INTRODUCTION 



In this paper we continue our study of the phase behavior of polydisperse fluid mixtures. In the previous studies 
IH S CJ El HI a scheme based on integral equation theory that allows one to include correlation between the particles 
on the level of the mean spherical approximation (MSA) has been developed. Usually the theoretical concepts used 
to study polydisperse systems treat such systems as a mixture with an infinite number of components, each of them 
characterized by a (continuous) variable x, which is distributed according to the distribution function f{x). In 
generalizing the formalism of the MSA for a mixture with a finite number of components to the polydisperse case, we 
r"| . obtained expressions for the structural and thermodynamic properties that can be expressed via a finite number of 
generalized moments of the distribution function f(x). Such models are called truncatable free energy (TFE) models 
and have the attractive feature that the formally infinitely many coexistance equations for phase equilibrium can be 
mapped onto a finite number of coupled, highly non-linear equations in these moments. Polydisperse mixture of hard- 
sphere Yukawa pHU and charged hard-sphere fluid [1 have been studied via this approach. In the case of Yukawa 
^ ■ mixtures application of the MSA is restricted to the systems with factorized version of the one- Yukawa interaction, i.e. 
the matrix of the coefficients describing the strength of the corresponding interaction is factorized into the product of 
two vectors. To remove this restriction recently high temperature approximation (HTA) for polydisperse multi- Yukawa 
hard-sphere fluid has been developed It was demonstrated that within HTA polydisperse multi- Yukawa hard- 



sphere model belongs to the family of truncatable free energy models. Unlike MSA, HTA includes correlations at the 
level of the reference system and in general it is less accurate than MSA, especially at lower densities However, 
HTA can be used to describe the phase behavior of a polydisperse hard-sphere Yukawa mixture with arbitrary number 
of non-factorized Yukawa tails. This feature makes the number of the systems amenable to HTA description much 
& \ larger than the number of the systems amenable to MSA description. 

In this paper we are focused on the phase behavior of polydisperse mixture of Yukawa hard-sphere chain molecules 
with polydispersity in the chain length. Our study is based on the extension of the recently proposed HTA for 
the model at hand. Laplace transforms of the site-site radial distribution functions (RDF) of the reference system, 
which are needed as an input, are calculated using product reactant Ornstein-Zernike approach (PROZA) 
appropriately reformulated in terms of the average total and direct correlation functions. For thermodynamical 
properties of the reference system we have used thermodynamic perturbation theory of Wertheim [T^ . Il3| , generalized 
. £^ | for polydisperse hard-sphere chain mixture. Similar as in the case of polydisperse hard-sphere Yukawa mixture 
polydisperse hard-sphere Yukawa chain mixture within present version of the HTA belongs to the family of TFE 
models. This feature allows us to use a scheme developed earlier []]■ 0. []| and to calculate complete phase diagram 
of the model (i.e., cloud and shadow curves as well as binodals) and chain-length distribution functions of the coexisting 
phases. 



II. HARD-SPHERE YUKAWA CHAIN MODEL 



We consider mixture of freely jointed tangent hard-sphere Yukawa chain molecules of N different chain species with 
densities p a and chain lengths m a . Here and below small Latin letters a, 6, c will denote species of chain molecules that 
take values 1..N, and small letters i,j,k, ... will indicate monomer species of the chain, which take values 1, ...,m Q . 
Monomers in each chain are assumed to be of the same size R a . In addition to the hard-sphere interaction Yukawa 
potential is acting between monomers, which are not connected directly. The resulting hard-sphere Yukawa potential 
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reads 

u y (r) = | _^ e _ z (o) (r _ il a i > )j ^ < ^ (1) 

where eo represents energy at contact, z(°) is Yukawa screening length and ro is distance unit. For simplicity reasons 
we have employed only one- Yukawa interaction, but our theory can be easily extended by multi- Yukawa potential, 
which can mimic large variety of realistic potentials. 

III. STRUCTURE PROPERTIES OF THE REFERENCE SYSTEM 

For the sake of simplicity instead of the regular site-site RDF gf b (r) we are using here average site-site RDF g ab (r), 
defined as follows: 



9 a \r) = -^— (2) 
m a mi, * — ' J 

ij 

These RDFs are calculated from the averaged version of the multidensity OZ equation, derived below. 



A. Averaged version of the multidensity OZ equation in the ideal chain approximation. 

Multidensity OZ equation for the multicomponent chain fluid in the ideal chain approximation is |qL llCt ITU]: 

<(r) = <» + W [c%{r') + A%(r')} p c a [h$(|r - r'|) + A$(|r - r'|)] dr' (3) 

where we have extracted from the correlation functions the terms with delta functions and denoted them by A^ b {r). 
The matrices used in are total correlation function matrix 

/ fc$oo(r) h$ 0A (r) hf b 0B (r) \ 
Kj(r) = h% A0 (r) h% AA (r) h? b AB (r) , (4) 
V h% m (r) h« b BA (r) hfl BB {r) J 



direct correlation function matrix 



c ij,Qo( r ) c ij.0A( r ) c ij.0B( r ) 

<>) = [ 4M r ) C ^A» S AB {r) ) , (5) 
,4,Bo( r ) c tsA(r) c% BB (r) 



density matrix 




(0) 



and the matrix containing delta-functions: 



\ / 



A$(r)=5 id+ i A« b (r) + | , (7) 

/ \ A ab (r) 



where 

A Qb fr) = 

4n Pa (R a 



A ab (r) = — j 5(r - R ab ), R ab = (R a + R b )/2. (8) 



Taking the average of both sides of equation Q we have 



m a m b 



m a m b 



E h #M = E c #( r ) + EE / [c^V) + A? fe c (r')]Pcak5(|r-r'|)+Ag.(|r-r'|)ldr'. (9) 



ij 13 ' ij ck 

By direct multiplication and summation, it can be shown, that 



£ £ A%(r')p c aA%(\r - r '|) = (m a - 2) £ A ac (r> c aA cfc (|r - r'|) 



where we have introduced new matrix 




A ab (r) = ( A ab (r) 
A ab (r) 



Similarly 



E A^ c (r')p c ahg(|r - r'|) = ]T(m c - l)m b A-(/)p c ah cfc (|r - r'| 



ij ck 



(10) 



(ii) 



(12) 



E E < fe c (r')Pc«Ag(|r - r'|) = £ m a (m c - l)c ac (r')p c «A cfc (|r - r'|), 

ij ck c 

where c cb (r) and h cb (r) are averaged direct and total correlation functions 



m a m b 



h ab (r) 



— E h ^ 



m a mb 



(13) 

(14) 
(15) 



It should be noted that A ab (r) has a Kroneker delta symbol 6 a b, thus sums over c in (|12() and 1|13|) can be removed. 
They are kept, because they will be used to assemble the OZ equation from the corresponding averaged quantities. 



h ab (r) = c ab (r)+J2Pcm c J 



c ac (r>h cb (|r _ r'|) + l -c ac (r')aA cb {\r - r'| 



^^A^(r')ah^.(|r - r'|) + I^A£(r>££(|r - r'|) 



Neglecting the terms with l/(m c ) 2 equation l|16ll can be factorized 



dr' 



(16) 



h ab (r)=c ab (r)+Y,Pcm c J 



c ac ( r ') + !^J-A ac (r') 



h cb (|r-r'|) + ^^A cb (|r-r'| 



dr', (17) 



and all the terms containing delta-functions can be included into corresponding correlation functions 



h ab (r) = c ab {r) + ^p c m c J c QC (r')ah cb (|r-r'|)dr', 



(18) 



where 
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c ab (r) = c ac (r) + — -A ac (r), (19) 

m c 

h ab (r) = b cb (r) + mc ~ 1 A cb (r), (20) 

m c 

are new averaged correlation functions which now depends only on chain species index, and are regarded to be the 
same for all beads of the chain. Also, new averaged contribution from sticky interaction was obtained, and it 
will be used now in formulation of PPY closure condition in terms of partial correlation functions h^(r) and c°^(r) 

htp(r) = -5 aO 6 0o , r<R ab , (21) 

c&(r) = 5ab{ma ~ 1] (8 aA S pB + 5 aB Sp A ) 6 (r R ab ) , r>R ab . (22) 
4np a m a (R ab ) 

As one can see the above closure conditions contain factor m ^ t ~ 1 which appears due to averaging of the OZ equation. 
It appears that the averaged version of the OZ equation (|18|l coincide with multidensity OZ equation for the fluid of 
associating particles forming a chains with the average length m a 0, 0] . 

B. Solution of the averaged OZ equation 

Solution of the set of the OZ equation 118fl with closure conditions 121fl and ljL'21) was performed employing Baxter 
factorization method ^t|- Factorizing this set of equations (|18f) we have 

-rc ab p (r) = [<$(r)] - 2^^p c m c ^(l - S l5 + ^ 8 so )- / ^ (r + t)dt, (23) 

r « S ac 



7,6 



for S ab <r< R ab , 



rh%{r)= [q%{r)]'- 27r$> c m c I> " &r« + M«>) T C 7 (*) (r - *) J$(|r - t\)dt, r>S ab , (24) 



7,(5 

where S^ 6 = ± (i?( a ) - and g is Baxter function 



9$( r ) = + ^) V + d%% (25) 



defined in interval S ab < r < R ab , with coefficients a£\b^ and A found to be 



1 1 i?( a )£ 2 



A 2 A 



2 



(26) 



4 A 2 

l(fi (a) )"6 1 R^R ab S ab & 
2 A + 4 A 2 



dSS = -o^-^ 1 + 7 > (28) 



(a) 1 (m a - 1) 

ay = t-(Sma + M> ( 29 ) 

b$ = -\^R (a \ (30) 

<o = -a&^S 06 , (31) 

df L = 0, (32) 

d ML = ~ A a ^ u(a) ( ~ ) {SmaSlB ~\~ ^LA^Ms) ; (33) 
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where & = 7r^ c p c m c (i? (c) )' and A = 1 - £3/6. 

The corresponding expressions for contact values of radial distribution function g^g(r) = <5 Q o^/?o + h^gir) are as 



follows: 



abl nab 

9oo\ H ~ 

„ab ( ryab 

9oM(R ab - 

9ML(R ab - 



1 1 R^R^b 



A 4 A 2 R ab ' 
1 i?^ (m a -\ 



AAR ab \ to, 
1 ijM (m h -l 



4Ai? Qb \ TO h 



to„ - 1 



(1 - ^ml) 



87rp a TO Q (R( a )) R ab \ m a 

Note that 5ml ( r ) contains also delta function 5 (r — R^j which describes neighbor hard spheres of the chain. 



(34) 
(35) 

(36) 

(37) 



C. Laplace transforms of the site-site radial distribution functions 

Integrating both sides of equation 1241) with J sab e~~ sr ...dr we have 
]T £ (S ca S Sa - 2n Pc m c [Q™ (a) + (1 - S SA ) Q a a c A (s) + (1 - S SB ) Q a a c B (s)]) G%{s) = 



(38) 



C 5 



/>OC 



fl ac / J gab £' 
7 5 + <5 7 o<5<5o) ^50^/30 / 9a 7 (*) — ) e ~ 



where Laplace transforms and their derivatives are denoted by capital letters: 



Q%{s) = / e- sr q%%{r)dr = e 



-sS ab 



S ab 







( 









(^ (b) ) + ^i(* (b) )) + + ^ (b) ) 



& (5) 



where 



(i-j a0 )(i-<W<W# (6) ) 

ab 

^ 2 (i? (b) ) 
<Mi? (b) ) 

^o(i? (b) ) 



+ ^»U (fl (6 ») 



5/50 (39) 



^0 (40) 



1 - + s 2 R^ '72 - e- 



1 - - e-^ ( ' 



1-e" 



(6) 



and 



= / e- sr g%%{r)dr. 



(41) 
(42) 
(43) 

(44) 



Similarly as in 

HI 13 HI, solution of equation (H reduces to inverting Jacobi like matrix. Note that here we 
have infinite number of vectors forming the matrix. To describe briefly our method of matrix inversion, let us start 
with equation (|38|l rewritten in following form: 



6 



(a) 



c 5 



1 RW 



2s 



2Trp b m b 



where Jacobi like matrix is given by 



^ - «L a) ^ c) - 4 a) 4 c) - E (4 a)(t) /j c)(t) + e )(t) ^ c)(t) ) • 



where 



4 C) = 2tx p c m c 



m(R {c) ) + ^vi(R (c) ) 



_ p -.flW /2 ja0g^ 

2A ' 



^=27rp c m c ^(^ c ))e sii(c, /2, 



5(o)W - p -^ (a) /2. TOa ~ 1 



/i c » t ^^(i?( c ))e^ ) /^ ct( l-^), 



'5 



;(c)(t) 



2m a i?( a ) 

<p (i? (c) )e sii(C)/2 ^t(l-^B), /or i = l..JV. 



The last four vectors e£ )[t) , /] c)(t) , g ( a a){t) and fc< c)(t) appear due to decomposition of the Kronccker delta 

Sab = E $a,Ab 



(46) 

(47) 
(48) 

(49) 

(50) 
(51) 

(52) 
(53) 

(54) 
(55) 



in d^p. This decomposition makes matrix Jacobi-like. In general it appears to be impossible to revert such 

matrix with arbitrary number of forming vectors e£°^*\ /j c ^*\ fja^^ and for t = 1..N without notion that 

some of these vectors are orthogonal and other can be composed back to form 5 ab : 



(efr>f<«>) = EE4 c)(p) 3 



<■><-/) ^ m P~ 1 ^ R (P)) S . 



C 7 



2m p i?(p) 



c 7 



c 7 



2m p R(pr UK ' pq 



(56) 
(57) 
(58) 
(59) 



The inverse matrix was searched in the following form: 

MX 1 = 6 xy 6 lX + a ab a^b{ y) + a cd c^d[ y) + a cb c^b^ + a ad a^d{ y) + ]T a™ (a%f[ vm + a a fh^ m ] 

i 

+ E i x) (°$ # )(t) + + E ( W> + °g# )(t) ) ^ + E KM* )(t) + °2# )W ) 4 W) 



EE 

tl *2 



(*l)(*2)-(x)(tl) f(l/)(*2) 



;(!/)(*.) + £ ag)^)^)^)^)^) + e E «il l)(i2) 4 a;)( * l) ^i y)( * 2) 



(60) 



tl t2 



tl *2 



+ EE4/ l)(t2) ^ x)(tl) ^ )( * 2 



ti t 2 
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where a a b, o.cd,...,or g j ' 2 are unknown coefficients that follow from the relation MM" 1 = I with I being unity matrix. 
We will skip these cumbersome calculations. 

After receiving unknown coefficients a the inverse matrix was multiplied by the right hand side of equation 145|) 
giving the following expressions for the Laplace transforms of the RDFs: 



^-sR a 



AD 



s 2 s ^ As 2 l-if( c ) (po(R^) 



4As 



+ - 



^ E w*»H W V0R (c) ) + E ^ mCyi 2A^ (il,(C)) (* W ~ R(C) ) ( Rib) RiC) 



+- 



2As 2 

c 

AD 1- K« 



(1 - <$ao) <5/30 



1 R 



• ab 



Ria) ( x ^ Ria) ^ 2 -Ew^ 2 e ( ^; 



?>o(*)) 



2s J 2As 2 



E 



4As 2 As 2 Lp {R {a) ) 
wp c m c ip 1 (R < -^) 



-i s Ew" C i? (c) V(i? ( 



2As 



ipo{R w ) J V 

1 i? Qb / 1 \ R^> 



+ - 



-sR" b K (b) 

AD l-ifW 



MR {b) ) 



2s 



2At 



R^R(»^ 2 1 R(c) 3., p(c) , 



+ E 



4As 

Trp c m c tpi 



+- 



AL> 1 - K( 

RWRW& 1 <p!(RW) (pi(R {b) ) 



4As 



^Tl - AW (! - *«o) (1 - W 

Ew™ci? (c) V(# (c) ) 



2 A ^ (i?W)^ (i? (b) ) 



2 A ^ (i? (a) ) ^o(i? w ) 



2A V" <*(J*0)'" , (#))jE^^^)^^« (c) ) 



+e~ 



2-Kp a m a 1- K( 



where 



TOq - 1 

2m a R(- a 



■<Po 



(i? (Q) ), 



^(i?^) - 1 + 

e(R (a) ) = i + 



2Jf( a ) 
1 - K^) ' 
2ir( a ) <^i(i? (a) ) 

i - P ip (R (a) Y 



(61) 

(62) 
(63) 
(64) 
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and 



D = 1 — 



;^ c m c i? (c) V(i? (c) ) 



^p t m^ 2 {R^)Q(R^)\l + 

t V c 

i n P tm t ip 2 (RW)R^(RM) j& - £ n Pc m c R^ 2 0(JJ< C >) J 



--^7rp c ,n c ^(i?( c ))ii?( c ) U(# (c) ) 



1 6^ (c 

2 4A 



9(i?( c )) 



(65) 



^( J RW)^ 7 rp c m c ^i(iiW)e(iiW) (R^ -R {t 



4A 

c y 

Presented above Laplace transform of partial RDF is reperesented as a sum of five terms that are multiplied by 
delta symbols and thus some of them vanish for a ceratin values of the monomer state indices a and [3. Laplace image 
of overall radial distribution function can be calculated by 



( ^ < overall( S ) ~ ^a^( S )- 

a/3 



(66) 



The first four terms of l|61(l contribute to intermolecular RDF, while the last term describes intramolecular RDF and 
contains Laplace transform of the delta-function. The intermolecular and intramolecular parts of the overall RDFs 
are related to partial RDFs by 



Gt ter {s) = G a b (s) + 2Gf (s)+2G a b A (s)+4G A b A (s), 

Gintra( s ) = 2 (G AB (s) — G AA (s)) . 



(67) 
(68) 



IV. THERMODYNAMIC PROPERTIES OF THE HARD-SPHERE YUKAWA CHAIN MIXTURE 

As in our previous paper Q , we use the HTA to obtain analytical expression for the Helmholtz free energy. Similarly 
as in the SAFT[13|,[21| we assume that free energy per unit volume / = A/V can be written as following 

/ = fid + fuS + fHTA + f chain, (69) 

where is the ideal gas contribution 

Pfid = ]T Pc m c [In (A 3 cPc ) - 1] , (70) 

c 

with A c being thermal de Broglie wavelength, fff s is hard sphere excess free energy and Jhta is key HTA term that 
accounts for attractive interaction between the chain monomers and contains Yukawa potential averaged with respect 
to the structure of the reference system. Finally the third term in l|69() contains contribution from the chain formation 



(71) 



where Y is hard-sphere cavity correlation function. 

The hard sphere excess free energy fff s and cavity correlation function Y are taken in Mansoori-Carnahan-Starling- 
Leland approximation!^, which for hard sphere system of equal sizes are given by 



0f. 



ex 
HS 



Y 



1 


l 2 R 3 




2ttA 


V 1+ 3A" 


1 


IR 3 


l 2 R e 


A + 


4A 2 


+ 72A3 ' 



(72) 
(73) 



where R is the hard sphere diameter, A = 1 — IR 3 /6 and / = ir ^ c p c m c is the moment of the chain number densities. 
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Key HTA term in l|69|l holds Yukawa attraction energy averaged with respect to the overall RDF of reference system 

2 , , poo 

ft] hta = 13- >~,7rp m„ y\ixp h mb / drr 2 uf{r)gl b verall (r). (74) 

n a b J ° 

Due to the choice of interaction potential JD integral in (|74|l can be substituted by overall RDFs Laplace transform 
of chain fluid, hence above formula can be written as follows 



PJhta = 



2f3e r e 



■S<»R 



Kp a m a *Pbm b e z(0) R G a b verall (z^), 



(75) 



where 



nab (U0)\ 
Lx overall \ z ) 



drre 



z (0 '>rab 



-dl( r )- 



(76) 



We will use obtained expression for the Laplace transform of the overall RDF (|61fl for HTA free energy contribution 
(|75|l . which now reads 



PI, 



HTA 



2f3e r 
irAD 



1 R R 3 <pAR) IR 4 
-v ■' 



R 3 



(I + 2v) 



+4v 2 



z (o) 2 Az(°) <p (R) 4Az(°) 2Az(°r 

R R '' 1 " N - • '■• " ,] 



(77) 



(l + 2v)R 3 tpijR) 

z (o) 2 ' 2 (o) tp (R) \ z (o) 2 ' 2z(°) J Az(°) M R Y ' 4Az(°) 1 2Az(°) ip (R) 



1 R 

+ 



R 6 ipAR) IR 
v 



MR) 
MR) 



R 



MR) l R 4 i + 2v({R 3 m r ) z(0)r2 \ (M r ) 



MR) 4Az(°) 



\MR) 



-R 2 MR) 



MR) 
MR) 



(l + 2v£$zW)R*MR) 



2Az(°) 



2(3e r Q 2v 
7T M R ) 



where 



MR) 



D = l- -MR) l + tv^-s 1 + — (I + 2v) + 



MR) 



R 3 



2A 



R 3 MR) 
A 2 



(l + 2v)2v^§ ]S 

MR) 



(78) 



(l + 2v) ( 1 + ^ ) + ( I + 2v 



2A 

I = 7T^^p C TO c , 



MR) 
'MR)' 



K, 



are the moments of the chain species distribution function, and 



(79) 
(80) 

, (81) 

(82) 
(83) 

To study the phase behavior, we will also need expressions for the pressure of the system and chemical potential of 
the chain molecules, which are derived from the above expression for Helmholtz free energy (|69|l - l|78(l 



MR {b) ) = 


l 


_ z (0) R (b) + 


z (0) 2 R(b) 2 / 2 






Z (0) 3 


MR [b) ) = 


l 


- zWRto - 






z (0) 2 




MR {b) ) = 


l 













PP = P^PcVc- Pf = PPHS + PPchain + PPHTA, 
c 

9W) 



Pp 



dp c 



(3p c ,id + PpThS + flpc,chain + (3p c ,HTA- 



(84) 
(85) 
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Here ^ C)j( j, fJ,^ HS , /J. c ,chain, Hc,hta and P c hain, Phta hold contributions from corresponding free energies f id , fff s , 
f chain and /htAi except for PPhs that holds sum of contributions from f ia - and fffs- These expressions for excess 
chemical potential and pressure can be formulated in terms of moments and their density derivatives 



where 



Placid 

0H 



ex 

c,HS 



PPhs 

PPchain 

PPhta 



m c \n{A 3 c p c } 
lR 3 l{c) 



7rA 



1 + 



3A 



+ 



1 



PVcchain = -(m c - l)ln{47r/9 c i? 2 r} - 

f3p c ,HTA : 



l 2 R e l(c) 
12ttA 2 

(I - irp) dY 



2 

3A 



dp c 



4/3e r e 



7T 



np c m c np a m a GZeraii ( 2(0) ) 



2pe r e z ' 0)R ^ ^ 9Gf veraU 

-^Tip a m a ^TTp b m b - 



(*(°)) 



7T 



Z 

Att 



LR3 
2A 



PR 6 
12A 2 



7rp 



7T 

2/3e r £ 



5Zi? 3 ; 2 i? 6 

+ 



12A 2 

R 



b 

l 3 R 9 
216A 2 

+ 1 



7T 



12A 3 
7rp a m a Y Kp b m b G\ 



ac 

overa 



h(* (0) ) 



2/3e r oe z(0)fi ^ ^ ^ 0G^ era 



,(*<°>) 



9p c 



ay 

dp c 

lc 



9A 2 

= 7rm c , 



144A 4 : 



5i? 3 / c . zi? 6 ; c . i 3 r 9 i c 

12A 3 

<9p c 
9w 



t'c = 



dp c 1 - A c 



(86) 
(87) 

(88) 
(89) 



(90) 
(91) 

(92) 



(93) 
(94) 
(95) 



and p = ^2 C p c is total number density. Because of the size expressions for p c ,HTA and Phta via the moments are 
shown in the Appendix. 

Thus, we can conclude, that our model belongs to the family of the TFE models, since its thermodynamic properties 
can be represented by finite number of (generalized) density distribution moments. 



V. TWO PHASE EQUILIBRIA 

When at constant T mother phase with the number density pi^ of the chain species c and volume separates 
into equilibrium daughter phases, then all the chemical potentials of every chain species and all pressures in each of 
the daughter phases are equal. Thus for the two-phase equilibrium we have 

^ ({p c 1} }) = Pf } ({ P ( c 2) }) , for every c (96) 

p(l) = p(2) ) 

where P&O and ({pi p) }) are pressure and chemical potential of phase p. We have explicitly specified the 

dependence on the species number densities |pc P ^| to indicate that they differ between phases. In addition particle 
number and volume conservation rules have to be satisfied 

/(°) = x WfP +x WfW (97) 
= x^vW+x^vW, (98) 
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together with normalization condition 

l = x^+x { - 2 \ (99) 

where x^ = /N^ for p = 1,2 is the chain molecules number fractions of phase p, = K^'p is the total 
chain molecule number in phase p, and /i^ = Pc/p^ — Nc /N^ is molecular species distribution function 

= 1. (100) 

c 

Assuming that we have the number of the chain species M, conditions (|96|) - (|100|l form closed set of M + 4 equations 
for M + A unknowns: fc , p^\ P^ 2 \ x^\ x^ 2 \ As in our previous work 6], to solve above equations, we will employ 
general scheme developed by Bellier-Castella et al 23] . The scheme is suitable only for truncated free energy models 
where very efficient transition to new unknown variables, represented by the moments, can be done. According to the 
scheme, the set of M + 4 equations, applied to our model, transforms to 4 equations in two unknown moments of the 
first phase i>W and two unknown densities of the daughter phases p^\ p^: 

p (1) = p (1) E /i 0) ^ ( T ; p (1 \i (1) v (1) ; p [2 \i (2) v^) , (ioi) 

c 

lW = pW £ l e fW He (T; p« ; p™ , Z< 2 V 2 )) , (102) 

C 

V W = p W J2 v c f®H c (T; pM,lMvW . ^2)^(2)^(2)^ (103) 

C 

P& (T;p<n,l<Uv<n) = P (2) {T;p^,l^v^) , (104) 

where 

/ \ (o^ - o^) A 12 (T- W iWrf 1 )' ( 2 ) /( 2 M 2 )) 
H c (T;pW,lWvU;p( 2 \lWvW) = -, [P x P > c [ '\ ' - ' P '- > ,(105) 

V ' 1 P^) + (pU - effi) A? (T lP WjWvU-,pP),lWv&i) 

Al 2 (T;pW,lW V W;pW,lWvW) = ^ exp [^ (2) (r; p (2) , Z (2 V 2) ) - pT W (t; P {1 \l (1) v^)] , (106) 

and pf c x = fi c — (iid,c is excess chemical potential. Note that two moments of the second phase and aren't 
unknown variables, because they can be found from the relation between moments of the phase 0, which are known 
in advance, and moments of the phase 1 

where m is one of the moments I or v. Set of equations presented above are used to locate the binodals. To obtain 
cloud and shadow curves, we consider the limit p^ — > p^°\ which employs that second (cloud) phase is very close to 
parent phase and thus there is only an infinitesimal amount of the first (shadow) phase . In such case it can be easily 
shown that 

Urn H c (T-p^\l^v^-p^ 2 \l^v^) = A? (t ; p« ,l^v^; p« ,l^v^) , (108) 
and therefore in equations i|101|) - (|lU3|) . H c should be replaced by A 12 . 



VI. RESULTS AND DISCUSSION 



We will present now phase diagrams for the two systems that differ from each other due to the mother distribution 
function /J . In contrast to the other studies of polydisperse fluids, we assume that distribution function depends 
on discrete species index c, rather than continuous one. But the species number should be large enough to expose 
properties of polydisperse systems. Each chain species c has length m c = c monomers, for c = 1..M. In most cases 
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FIG. 1: Phase diagram (T* vs p*) of polydisperse chain fluid with D m — 0.5 beta distribution function. Cloud curve, shadow 
curve are represented by solid and dashed curves respectively. Dotted curve denotes critical binodal, and dash-dotted line 
represents monodisperse binodal. Also five points on shadow curve Ai,Bi,Ci,Di,Ei together with corresponding equilibrium 
points A2,B2,C2,D2,E2 on cloud curve, mark points for which distribution functions are presented on [5] 
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FIG. 2: The same as figure but for Shultz distribution with z = 1. Corresponding distributions for specified points are 
presented on llOl 

we consider the systems with M = 100. For the analysis of the cutoff effects we present results for 200 and 300 chain 
species. 

In the following we will consider two systems with two different distribution functions. The first mother distribution 
function was chosen to be Beta-distribution 

^-ssflterVs)" (io9) 

for c = 1..M, where B(a,0) is beta function and a and (3 are related to mean chain length (mo = (m) ) and mean 
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FIG. 3: Phase diagram (T* vs 77 = (tt/6)R 3 ^2 c m c p c ) of polydisperse chain fluid with D m — 0.5 beta distribution function. 
The curves denote the same as inQand|3 



square length (m 2 ) by 



a 



D r , 



M - m (1 + D m ) 



AID 
M - m 

too 
(to 2 ) 
toq 2 



a, 



mo 



(110) 
(111) 

(112) 



with a m = Y \ ( m — *! m )) / being standard chain length deviation. 

The second distribution function was chosen to be represented by the Shultz distribution 

2+1 

c z exp 



f (o) = z. 

J C 1 



1 fz+l 



toq 



z + 1 



toq 



(113) 



for c = 1..M, where z is dispersity controlling parameter. 

In both cases the dispersity of the distribution functions were chosen to be equal, with standard deviation a m = 11.31 
that corresponds to D m = 0.5 for beta distribution and z = 1 for Shultz distribution, and with average chain length 
too = 16 for both distributions. 

The corresponding phase diagrams for both systems are demonstrated on figures ^ and [21 in terms of the number 
densities and on figures |3] and 0] in terms the volume fraction rj. Phase diagrams on figures |3] and 0] are almost 
indistinguishable, but figures HO have one significant difference: liquid part of the shadow curve for the system with 
Shultz distribution has smaller density than that of them phase diagram with Beta distribution. This difference of 
the phase behavior can be explained by inspecting distribution functions on figures 151 and HU1 for corresponding points 
marked on the phase diagrams on figures ^ and [21 Liquid shadow curves for both systems consist of the chains longer 
than average, in contrast to the gas shadow curves, that are mainly composed of shorter chains. As one can see for 
the system with Shultz distribution the chains in liquid phase are on average larger than for the systems with Beta 
distribution. The origin of these effects lies in the difference of the tails for the mother phase distributions (either 
Beta or Shultz). For Shultz distribution the tail is "fatter" than for Beta distribution. Being more precise at cutoff 
distance m cuto ff — M = 100 the Shultz distribution tends to finite value as opposed to Beta whose value tends to 0. 
These differences in compositions of the two liquid shadow phases are the reason for the discrepancy between phase 
diagrams. Note, that volume fractions are close, only average chain lengths and as a result number densities differ. 

Figures O and El show shadow curves and critical binodals in terms of the average chain length (to). Average chain 
length of the cloud curve is always 16. In both cases liquid part of the shadow curve moves to longer chains with the 
temperature decrease, but for Shultz distributed system this move is more rapid. 
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FIG. 4: Same as [31 but for Shultz distribution with 8 = 1. 



Figures [7| and [H] expose shadow curves and critical binodals in terms of the chain length standard deviation a m . 
Similarly to the average lengths liquid part of the shadow curve is on the right and as the temperature decreases it 
moves to larger chain length deviations. However after reaching some maximum it starts to decrease. This decrease 
corresponds to the formation and growing of the peak on the distribution function in the region of long chains. Beta 
distributed system shows noticeable peak for the point E\ on the phase diagram, its distribution function can be 
found on figure With the temperature decrease the peak grows and approaches the right margin M = 100 of chain 
lengths. The peak for the shadow phase of Shultz distributed system has irregular form, its maximum is always at 
c = M = 100 (see figure^Jfor point J x ) and it also grows with the temperature decrease. These effects indicate that 
for the low temperatures liquid shadow phase becoming more and more monodisperse and is composed mainly of long 
chains. The same happends for the gas shadow phase, but it "likes" more short chains than long. 

To study conditions at which one can observe irregular peak formation, we examined six additional phase diagrams 
of Shultz distributed system with z — 5, z — 10 and at three values of the cutoff chain lengths, i.e. the number of 
species m cu toff = M = 100, 200, 300. They can be found at figure 11. We can deduce that increase of the cutoff chain 
length lowers temperature of the peak formation. However this lowering is rather negligible when compared to the 
change of the peak appearance temperature due to the change of the mother distribution polydispersity. Thus, to get 
the peak on the distribution function of the liquid shadow phase at high enough temperature one should study the 
system with high enough polydispersity. 

We also examined the change of the critical point position due to the change of the mother phase distribution 
function polydispersity ffigures IT21 and IT3|) . For both systems critical point moves to the higher temperatures and to 
the lower volume fractions with the increase of the distribution function width. 

In our previous work |6|], which is dedicated to the Yukawa hard spheres with size polydispersity we concluded 
that hard sphere repulsive interactions plays a major role in controlling fractionation effects. In the present case the 
situation is different: strong intramolecular attraction may change effective radius of the chain molecules, which are 
flexible and this may have influence on the phase behavior. But our theory takes structural properties of the system 
only on the level of hard sphere chains without attraction, and hence it can't capture the effects of effective radius 
change. 

To verify the accuracy of our theory we need to perform molecular MC simulations similar to those made by Wilding 
atalHHii. 



VII. CONCLUDING REMARKS 



We have shown that polydisperse Yukawa chain mixture with chain length polidispersity treated within HTA 
belongs to the TFE class of the models, i.e. its thermodynamical properties were shown to depend on the finite 
set of generalized distribution function moments. Exploiting this property we represent coexistence relations as a 
set of four nonlinear equations for unknown moments of daughter distribution functions. To illustrate how tails of 
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FIG. 5: Shows shadow and critical binodal curves for system with Beta distribution, denoted respectively by dashed and dotted 
lines, in terms if T* and phase average chain length (m). 




90 100 



FIG. 6: Same as figureQ3 but for Shultz distribution. 



the distribution functions can impact the phase behavior these equations were solved and the phase diagrams of the 
systems with two different mother phase distribution functions with the same dispersion characteristics (mean chain 
length and standard deviation) and different asymptotics of the distribution function tails were studied. The obtained 
results show strong dependence of the shadow curve position on the asymptotic of the mother phase distribution 
function. 

We also demonstrated that temperature of the appearance of the peak on the distribution function of the shadow 
phase strongly depends on the system polydispersity. 



APPENDIX A: EXPRESSIONS FOR THE CHEMICAL POTENTIAL AND PRESSURE 



The quantities given below can be simplified, but they are kept her in the form that resembles the structure of the 
overall RDFs Laplace transform. 



16 




FIG. 7: The same as [3] but reexpressed in terms of chain length standard deviation a m = \J (m — (m) ) . 
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FIG. 8: Same as figure|7| but for Shultz distribution. 
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FIG. 9: Distribution functions f c = p c /p of Beta distributed system of five points, marked on0 along shadow curve, together 
with cloud curve (mother) distribution. 
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FIG. 10: Same as figure^l but for system with Shultz distribution, and phase diagram, containing the points, given on figure|21 
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FIG. 11: Figure shows shadow curves of system with Shultz distribution for z = l,z = 5 and z = 10 and different cutoffs of 
distribution function, i.e. different maximal species chain lengths 100,200 and 300. The figure is displayed in terms of T* and 

standard deviation of chain length a m = \J(m— (m) ) 2 . The decrease of chain length deviation with decrease of temperature 
that is located below chain length deviation maximum, corresponds to formation of distribution function peak near the cutoff 
distance. The peak appears in liquid branch of shadow curve. 



3.1 



3.05 - 



2.95 



2.85 



2.8 



T* 



t7m fa m 3 i 3 22.63 

■ a m = 21.47 

■ cr m = 20.24 

■ a m = 18.93 



I <r m = 17.53 

■ o m = 16.00 

■ a m = 14.31 



I cr m = 12.39 
■ a m = 11.31 

■ a m = 10.12 



I a m = 7.16 
■ o m = 5.06 



T] 



0.086 0.088 



0.09 0.092 



0.094 



0.096 0.098 



0.1 



FIG. 12: Critical temperatures and volume fractions rj for variety of Beta distributions dispersions. 
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FIG. 13: Critical temperatures and volume fractions r\ for variety of Shultz distributions dispersions. 



Contribution from HTA free energy term to chemical potential reads 
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is contribution from intramolecular Yukawa interactions to free energy (|78|l and dD/dp c has the following form 



dD 

dp c 



3A 2 



I + 2v^Mz^ 



MR) 

2v 7— rS 



( 



MR) 

R 6 L 



MR) ) V12A 



(I + 2v) 



' 2A 
2A 



[(J c + 2v c ) 2v + (l + 2v) 2v r ] 



(I + 2v) 
{lc + 2v ( 

(0) 



-MR) 



MR) 



2A 



{I + 2v) 



3A 3 



R 3 



-^MR)R 



IR 3 

2^M 1+ 2A 



MR)' 

(l + 2v) 



3A : 



l MR)R 



MR)' 
m 3 



R 3 L IR% 



2A 12A 2 



(1 + 2v) ( 1 + + + 2v 



. MR) 
MR)' 



(A3) 



(0) 



MR) 



Corresponding pressure expression is given by 
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where ^ c p c §^- can be obtained from IA3I by simple replacing of all l c and v c by I — ^ c p c l c and t> = ^ c p c v c 
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